#include "ImageGray.h"
#include "tp3.h"
#include "fft2d_util.h"
#include "ipFFT2D.h"
#include "DiskXFile.h"
#include "xtprobeImageDiskXFile.h"


fft2dTf::fft2dTf(char *name) {
	spectrum = NULL;
	initial_image.reset(new Torch::ImageGray());
  	//ImageDiskXFile *image_file = new xtprobeImageDiskXFile(name,"r");
	std::tr1::shared_ptr<Torch::ImageDiskXFile> image_file(
			new Torch::xtprobeImageDiskXFile(name, "r"));


	file_name = (char*) malloc (sizeof(char)*longueur(name)+1);
	
	strcpy(file_name,name);

  	//chargement de l'image
	initial_image->loadImageXFile(image_file.get());
}

void fft2dTf::fft() {
	
	int new_width = (int) pow(2.0, ceil(log2((double)initial_image->width)));
  	int new_height = (int) pow(2.0, ceil(log2((double)initial_image->height)));

	power_2_image.reset(new Torch::ImageGray(new_width, new_height));

  	pasteGray(initial_image->data, 0, 0, initial_image->width, initial_image->height, power_2_image->data, new_width, new_height);
  	power_2_image->updatePixmapFromData();	

	// calculer fft
	fft2d.reset(new Torch::ipFFT2D(power_2_image->width,power_2_image->height));
	fft2d ->process(power_2_image.get());

	fft_real = fft2d ->seq_out->frames[0];
	fft_imaginary = fft2d->seq_out->frames[1];

	spectrum = new float[fft2d->seq_out->frame_size];

	return;
}

void fft2dTf::reverseFft(char *name) {

	//char reverseImageName[] = {'_','i','n','v','e','r','s','e','.','p','g','m','\0'};
	
	// calculer fft inverse
	ifft2d.reset(new Torch::ipFFT2D(power_2_image->width,
			power_2_image->height, true));
	ifft2d->process(fft2d->seq_out);

	// sauvegarde de l'image reconstruite
	ImageGray *image_out = new ImageGray(initial_image->width, initial_image->height);

cropGray(ifft2d->seq_out->frames[0],power_2_image->width,power_2_image->height,0,0,initial_image->width, initial_image->height, image_out->data);
 	image_out->updatePixmapFromData();
	
	image_out->save(newName(file_name,name));
	return;
}

//
// fromTab1D_to_Immage2d
//
// Cette fonction transforme un tableau 1D en une image 2D et lùenregistre sur le disque
//
// Entrees
// 
// improve : rehausse le contenu de l'image par une transformation logarithmique
// swap : echange les quadrants pour mettre l'origine au centre de l'image
//choice de l'operation a realiser 
//
// Sortie
// une image ImageGray de taille width x height
//




void fft2dTf::fromTab1D_to_Immage2d(bool improve, bool swap, int choice)
{
	if(choice == 1)
	{
		char realName[] = {'_','f','r','e','a','l','.','p','g','m','\0'};
		ImageGray *imgreelle = (ImageGray*) Tab1D_to_Image2D (fft_real, power_2_image->width, power_2_image->height, improve,swap);
		DiskXFile *image_file1 = new DiskXFile(newName(file_name,realName),"w");
    		imgreelle->saveXFile(image_file1);
		delete image_file1;
	}
	else if(choice == 2)
	{
		char imaginaryName[] = {'_','f','i','m','a','g','i','n','a','r','y','.','p','g','m','\0'};
		ImageGray *imgimaginaire = (ImageGray*) Tab1D_to_Image2D (fft_imaginary, power_2_image->width, power_2_image->height, improve,swap);
  		DiskXFile *image_file2 = new DiskXFile(newName(file_name,imaginaryName),"w");
  		imgimaginaire->saveXFile(image_file2);	
		delete image_file2;
	}
	else if(choice == 3)
	{	
		char spectrumName[] = {'_','f','s','p','e','c','t','r','e','_','f','.','p','g','m','\0'};
		compute_spectre();
		ImageGray *imgspectre = (ImageGray*) Tab1D_to_Image2D (spectrum, power_2_image->width, power_2_image->height, improve,swap);
  		DiskXFile *image_file3 = new DiskXFile(newName(file_name,spectrumName),"w");
  		imgspectre->saveXFile(image_file3);
		delete image_file3;	
	}
	else // on transforme tous les tableaux de reels en image 2D 
	{// il s'agit des tableaux  spectrum, fft_imaginary et fft_real 

		char realName[] = {'_','r','e','a','l','.','p','g','m','\0'};
		ImageGray *imgreelle = (ImageGray*) Tab1D_to_Image2D (fft_imaginary, power_2_image->width, power_2_image->height, improve,swap);
		DiskXFile *image_file1 = new DiskXFile(newName(file_name,realName),"w");
    		imgreelle->saveXFile(image_file1);
		
		char imaginaryName[] = {'_','i','m','a','g','i','n','a','r','y','.','p','g','m','\0'};
		compute_spectre();
		ImageGray *imgimaginaire = (ImageGray*) Tab1D_to_Image2D (fft_real, power_2_image->width, power_2_image->height, improve,swap);
  		image_file1 = new DiskXFile(newName(file_name,imaginaryName),"w");
  		imgimaginaire->saveXFile(image_file1);
		
		char spectrumName[] = {'_','s','p','e','c','t','r','e','_','f','.','p','g','m','\0'};
		
		compute_spectre();
		ImageGray *imgspectre = (ImageGray*) Tab1D_to_Image2D (spectrum, power_2_image->width, power_2_image->height,improve,swap);
  		image_file1 = new DiskXFile(newName(file_name,spectrumName),"w");
  		imgspectre->saveXFile(image_file1);
		delete image_file1;			
	}
	
}


void fft2dTf::applayfiltre1(float factor)
{
	real *filtre = (real*) DefinirFiltre1 (power_2_image->width, power_2_image->height,factor);
	int taille = power_2_image->width*power_2_image->height;
	AppliquerFiltre1 (fft2d->seq_out->frames[0],filtre, taille);
	AppliquerFiltre1 (fft2d->seq_out->frames[1],filtre, taille);

	return;
}

void fft2dTf::compute_spectre()
{
	for(int k=0;k<fft2d->seq_out->frame_size;k++)
	spectrum[k] = sqrt(fft2d->seq_out->frames[0][k] * fft2d->seq_out->frames[0][k] + fft2d->seq_out->frames[1][k]* fft2d->seq_out->frames[1][k]);	
}


void fft2dTf::applayfiltre2(float factor, bool inverse)
{
	real *filtre = (real*) DefinirFiltre2 (power_2_image->width, power_2_image->height,factor,inverse);
	AppliquerFiltre2 (fft2d->seq_out->frames[0],filtre, power_2_image->width*power_2_image->height);
	AppliquerFiltre2 (fft2d->seq_out->frames[1],filtre, power_2_image->width*power_2_image->height);
	return;
}

/*
La fonction ci-dessous est utilisée pour corriger les effets du filtre 1 sur une image,
Le principe est le suivant : on recherche les points correspondant au sinus dans le spectre de fourier
la valeur des pixels a ces point est mise à 0, puisaue nous voulons autant que possible retouver exactement 
l'image de départ, nous mettons aussi les points aux voisins de celui représentant le sinus à 0 
*/
void fft2dTf::removeNoise(float facteur){

	int facty = static_cast<int> (static_cast<float> (power_2_image->height/ 2) * facteur);
	int factx = static_cast<int> (static_cast<float> (power_2_image->width/ 2) * facteur);
	
	real *filter;
	real *filter1;
	
	//on met tous les pixels au environs du point representant le sinus a 0
	if(facteur <= 0.5 ){
		filter = (real*) DefinirFiltre3 (power_2_image->width, power_2_image->height,facteur, false,(facty * power_2_image->width),factx);
		filter1 = (real*) DefinirFiltre3 (power_2_image->width, power_2_image->height,facteur, false,((power_2_image->height - facty) * power_2_image->width),(power_2_image->width - factx)); 
	}
				
	if(facteur > 0.5 ){
		filter = (real*) DefinirFiltre3 (power_2_image->width, power_2_image->height,(1+ facteur)/2, false,(facty * power_2_image->width),factx);
		filter1 = (real*) DefinirFiltre3 (power_2_image->width, power_2_image->height,(1+ facteur)/2, false,((power_2_image->height - facty) * power_2_image->width),(power_2_image->width - factx)); 
	}
	
	AppliquerFiltre3 (fft2d->seq_out->frames[1], filter,power_2_image->width* power_2_image->height);
	AppliquerFiltre3 (fft2d->seq_out->frames[0], filter,power_2_image->width* power_2_image->height);
	
	AppliquerFiltre3 (fft2d->seq_out->frames[1], filter1,power_2_image->width* power_2_image->height);
	AppliquerFiltre3 (fft2d->seq_out->frames[0], filter1,power_2_image->width* power_2_image->height);	
	return;
}

char* fft2dTf::newName(char*oldname,char*additionName)
{
	//on soustrait 4 pour fairepartir les informations du type de l'image
	int l = longueur(oldname) - 4;
	int size = l + longueur(additionName) + 1; 

	char* nouveau = (char*)malloc(sizeof(char)*(l+longueur(additionName) + 1));
	int i = 0;
	for(i = 0;i<l;i++)
		nouveau[i] = oldname[i]; 
	//printf("\t%s\n",nouveau);
	for(int t=0,k=i; k< size;t++,k++)
		nouveau[k]=additionName[t];

	return (nouveau);
	
}


int fft2dTf::longueur(char* chaine)
{

	int k =0;
	for(k=0;chaine[k]!='\0';k++);
	
	return(k);
		
}
// DefinirFiltre3
// cette fonction definie un filtre permettant de reduire le briut sinus dans l'image
//les elements xcord et ycord sont les coordonnees du point autour duquel on veut appliquer le filtre  
real* fft2dTf::DefinirFiltre3 (int width, int height, float facteur, bool inverse,int xcord,int ycord)
{
	real val1 = 0.;
	real val2 = 1.;
	int size = width * height;
	real *filtre = new real [size];
	for (int i=0; i<size; i++) filtre[i] = val2;

	int height2 = height/2;
	int width2  = width/2;

	int radius = int(float(height2)*facteur);
	if (width2 < height2)
		radius = int(float(width2)*facteur);

	int radius2 = radius*radius;

	for (int y=0; y<=radius; y++)
		for (int x=0; x<=radius; x++)
		{
			if (((y*y) + (x*x)) <= radius2)
			{
				filtre[(y*width + ycord) + x +xcord] = val1;
				filtre[(y*width + ycord) - x + xcord] = val1;
				filtre[(xcord + x) -(y*width - ycord)] = val1;
				filtre[(xcord - x) -(y*width - ycord)] = val1;
			}
		}

   return filtre;

}
//Applique le filtre fourni par la fonction DefinirFiltre3

void fft2dTf::AppliquerFiltre3 (real *fft, real *filtre2, int size)
{
	for (int i=0; i<size; i++)
		if (filtre2[i] < 0.1) fft[i] = 0;

	return;
}


fft2dTf::~fft2dTf() {
	if (spectrum != NULL) {
		delete[] spectrum;
	}
}
